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We determine how the differences in the treatment of the sub-filter-scale physics affect the prop- 
erties of the flow for three closely related regularizations of Navier-Stokes. The consequences on 
the applicability of the regularizations as sub-grid-scale (SGS) models are also shown by examining 
their effects on super-filter-scale properties. Numerical solutions of the Clark— a model are com- 
pared to two previously employed regularizations, the Lagrangian- Averaged Navier-Stokes a— model 
(LANS— a) and Leray— a albeit at significantly higher Reynolds number than previous studies, 
namely Re ~ 3300, Taylor Reynolds number of Rx ^ 790, and to a direct numerical simulation 
(DNS) of the Navier-Stokes equations. We derive the Karman-Howarth equation for both the 
Clark— a and Leray— a models. We confirm one of two possible scalings resulting from this equation 
for Clark— a as well as its associated energy spectrum. At sub-filter scales, Clark— a pos- 
sesses similar total dissipation and characteristic time to reach a statistical turbulent steady-state 
as Navier-Stokes, but exhibits greater intermittency. As a SGS model, Clark— a reproduces the 
large-scale energy spectrum and intermittency properties of the DNS. For the Leray— a model, in- 
creasing the filter width, a, decreases the nonlinearity and, hence, the effective Reynolds number is 
substantially decreased. Therefore even for the smallest value of a studied Leray— a was inadequate 
as a SGS model. The LANS— a energy spectrum ~ , consistent with its so-called "rigid bodies," 
precludes a reproduction of the large-scale energy spectrum of the DNS at high Re while achieving 
a large reduction in numerical resolution. We find, however, that this same feature reduces its 
intermittency compared to Clark— a (which shares a similar Karman-Howarth equation). Clark— a 
is found to be the best approximation for reproducing the total dissipation rate and the energy 
spectrum at scales larger than a, whereas high-order intermittency properties for larger values of a 
are best reproduced by LANS— a. 

PACS numbers: 47.27.ep; 47.27.E-; 47.27.Jv; 47.50.-d 

I. INTRODUCTION 

Nonlincaritics prevail in fluid dynamics when the Reynolds number, Re, is large p7j . For geophysical flows, the 
Reynolds number is often larger than 10^ and for some astrophysical flows values of Re « 10^^ is not unreasonable. 
The number of degrees of freedom (dof) in the flow increases as Re^^^ for i?e > 1 in the Kolmo gorov framework 
jsil . [3^ . [33! (hereafter K41). Such a huge number of dof makes direct numerical simulations (DNS) of turbulence 
at high Re infeasible on any existing or projected computer for decades to come. Because of this intractability, 
simulations of turbulence are always carried out in regions of parameter space far from the observed values, either 
with: (a) an unphysical lack of scale separation between the energy-containing, inertial, and dissipative ranges while 
parameterizing the missing physics, or (b) a study of the processes at much smaller length scales, often with periodic 
boundaries (unphysical at large scales but used under the hypothesis of homogeneity of turbulent flows). Clearly, 
modeling of unresolved small scales is necessary. 

Given the nonlinear nature of turbulent flows and the ensuing multi-scale interactions, the physics of the unresolvable 
scales may not be separable from the properties (e.g., statistics) of the resolvable large scales. However, two main 
approaches have been developed over the years to model the effects of the unresolvable small scales in turbulence on 
the scales resolved in the simulations. The first approach is Large Eddy Simulations (LES, see [1^). LES is widely 
used in engineering, in atmospheric sciences, and to a lesser extent in astrophysics. However, in the LES approach, 
the Reynolds number is not known. Instead, one attempts modeling the behavior of the flow in the limit of very large 
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Re. As the Kolmogorov assumption of self-similarity is known to be violated (e.g., by intermittency (sol. [soj and by 
spectral non- locality [l[), the value of Re can play an important role, e.g., in the competition between two or more 
instabilities [4^ . Therefore, another approach models the effects of turbulence at higher Reynolds numbers than are 
possible with a DNS on a given grid, by using a variety of techniques that can be viewed as filtering of the small scales 
(the so-called sub-grid-scale (SGS) models). 

A rather novel approach to modeling of turbulent flows employs regularization modeling as a SGS model [ill, fl^. fiol. 
HOj HH, . Unlike closures which employ eddy- viscosity concepts (modifying the dissipative processes), the approach 
of regularization modeling modifies the spectral distribution of energy. For this reason, they retain a well-defined 
Reynolds number. Existence and uniqueness of smooth solutions can be rigorously proven, unlike many LES models 
(e.g., eddy- viscosity) , as well as the fact that the subgrid model recovers the Navier-Stokes equations in the limit 
of the filter width going to zero. Their robust analytical properties ensure computability of solutions. These same 
properties reopen theoretical possibilities first explored by Leray when he proved the existence (but not smoothness, 
or uniqueness) of solutions to the Navier-Stokes equations in M", n=2,3 using the Leray model [sj- This treatment 
of the small scales, then, enforces a precise type of regularization of the entire solution which may be studied as an 
independent scientific questi on ( as compared to either LES or SGS modeling). 

Geurts and Holm ^ HO, [SI be gan using the Leray model with (a three-point invertible approximation of) an 
inverse-Hclmholtz-operator filter of width a. Later it was dubbed Leray— a and an upper bound for the dimension 
of the global attractor was established [l3|- The global existence and uniqueness of strong solutions for the Leray 
model is a classical result [IJ- Leray— a has been compared to DNS simulations on a grid of = 192^ in a doubly- 
periodic compressible channel flow domain [l^, [2^, [2li • Its performance was found to be superior to a dynamic mixed 
(similarity plus eddy- viscosity) model (with an even greater reduction in computational cost). However, it possessed 
a systematic error of a slight over-prediction of the large scales accompanied by a slight under-prediction of the small 
scales. It possessed both forward- and back-scatter, but exhibited too little dissipation. 

The Leonard tensor-diffusivity model [s^ (sometimes known as the Clark model [l^) is the first term of the 
reconstruction series for the turbulent sub-filter stress for all symmetric filters that possess a finite nonzero second 
moment. This leading order approximation of the subgrid stress is thus generic [3, [12 ■ In a priori testing, it 
reconstructs a significant fraction (> 90%, but not all) of the subgrid stress, provides for local backscatter along 
the stretching directions while remaining globally dissipative, and possesses a better reconstruction of the subgrid 
stress than the scale-similarity model. Used as a LES, in a posteriori testing the Leonard tensor-diffusivity model 
required additional dissipation (a dynamic Smagorinsky term) to achieve reasonable gains in computation speed for 
3D periodic flows and for channel flows [m|. The Leonard tensor-diffusivity model does not conserve energy in the 
non-viscous limit. Cao, Holm, and Titi [6| developed a related (conservative) subgrid model which they dubbed 
Clark— a. The Clark— a model applies an additional inverse-Helmholtz filter operation to the Reynolds stress tensor 
of the Clark model. The global well-posedness of the Clark— a model and the existence and uniqueness of its solutions 
were demonstrated, and upper bounds for the Hausdorff (d/f ) and fractal (cJf) dimensions of the global attractor were 
found This model has yet to be evaluated numerically. 

The third regularization model we will consider is the incompressible Lagrangian-averaged Navier-Stokes (LANS— a, 
a— model, also known as the viscous Camassa-Holm equation H, [^, [13, [O, [1^ ) ■ It can be derived, for instance, by 
applying temporal averaging to Hamilton's principle, where Taylor's frozen-in turbulence hypothesis (the only approx- 
imation in the derivation) is applied as the closure for the Eulcrian fluctuation velocity in the Reynold decomposition, 
at linear order in the generalized Lagrangian mean description [2^ . [26l . [27| . In this derivation, the momentum- 
conservation structure of the equations is retained. For scales smaller than the filter width, LANS— a reduces the 
steepness in gradients of the Lagrangian mean velocity and thereby limits how thin the vortex tubes may become 
as they are transported, while the effect on larger length scales is negligible [ll|. LANS— a may also be derived by 
smoothing the transport velocity of a material loop in Kelvin's circulation theorem [l6| . Consequently, there is no 
attenuation of resolved circulation, which is important for many engineering and geophysical flows where accurate pre- 
diction of circulation is highly desirable. An alternative interpretation of the a— model is that it neglects fluctuations 
in the smoothed velocity field, while preserving them in the source term, the vorticity [31 . 

LANS— a has previously been compared to direct numerical simulations (DNS) of the Navier-Stokes equations at 
modest Taylor Reynolds numbers (i?A « 72 [H], i?A « 130 [IH, and Rx ^ 300 Q^). LANS —a was compared to a 
dynamic eddy- viscosity LES in 3D isotropic turbulence under two different forcing functions (for Rx « 80 and 115) 
and for decaying turbulence with initial conditions peaked at a low wavenumber (with Rx « 70) as well as at a 
moderate wavenumber (with Rx « 220) LANS— a was preferable in these comparisons because it demonstrated 
correct alignment between eigenvectors of the subgrid stress tensor and the eigenvectors of the resolved stress tensor 
and vorticity vector. The LES effectiveness of the LANS— a and the Leray— a regularization models relative to 
eddy-viscosity and the dynamic mixed model (similarity plus eddy-viscosity) have already been demonstrated in a 
turbulent mixing shear layer (with Re « 50) [3, [2lj. LANS— a was found to be the most accurate of these three 
LES candidates at proper subgrid resolution, but the effects of numerical contamination can be strong enough to 
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lose most of this potential. While LANS— a has the greatest grid-independent accuracy of the three models, it also 
requires the greatest resolution. From the LES perspective, this could pose some limitations on the practical use 
and application of LANS— a for high Re cases. Indeed, recent high-resolution simulations of LANS— a showed that 
energy artificially accumulates in the sub-filter-scales, giving as a result only a modest computational gain at very 
high Reynolds number [i^ . 

We propose to pursue these previous studies of Leray— a and LANS— a further at higher Reynolds number, and 
to use them as a benchmark for evaluation of Clark— a. One goal is to contrast the sub-filter-scale physics of the 
three models to determine the relevant features from which to build improved models. As the three regularizations 
are related via truncation of sub-filter stresses, such a comparison can be illuminating. For LANS— a, the predicted 
sub-filter-scale spectra is ~ [l^. This scaling has been observed to be subdominant to an energy spectrum ^ 
which corresponds to "enslaved rigid body" or "polymerized" portions of the fiuid [i^ . The sub-filter scaling observed 
in the third-order structure function corresponded to the predicted ^ fc^^ scaling of the energy spectrum. However, 
regions were observed in the flow where no stretching was acting in the sub-filter scales. These regions, which give no 
contribution to the energy cascade, and hence do not affect the third order structure functions, are responsible for the 
scaling in the LANS— a energy spectrum. For Clark— a, the correct time scale for vortex stretching is difficult 
to determine and its spectrum is found to range between and ~ @. Leray— a has the same difficulty 

and the spectrum can range between ^ k^^^^ and ^ k~^^^ [l^. The determination of these scaling laws is needed 
to quantify the computational gain if each model is to be used as a SGS model. As a result, we seek to determine 
empirically the sub-filter scale spectra. Our second goal is to evaluate the applicability of these three regularizations 
as SGS models. This is accomplished both through prediction of computational gains from observed sub-filter-scale 
properties and through directly testing their capability to predict super-filter-scale properties at high Re. 

We present the three models and describe how they are related, derive the Karman-Howarth equation for Clark— a 
and Leray— a (from which exact scaling laws for third order quantities follow), and review theoretical predictions of 
incrtial range scaling in Section |TT1 We examine the sub-filter-scale properties of the three regularizations in Section 
mil We first compute a fully resolved DNS of the Navier-Stokes equations at a resolution of 1024^ (i/ = 3 • 10^^, 
Re « 3300, and R\ w 790). We then perform model runs with the exact same conditions at a resolution of 384^. 
We take a to be l/13th the box size, which was found in an earlier study to be large enough to exhibit both Navier- 
Stokes and sub-filter-scale LANS— a dynamics [i^. This large filter case is important because it gives insight into the 
behavior of the models at scales much smaller than the filter width without requiring higher resolution than is feasible. 
We compare the three regularizations as subgrid models in Section HVl Guided by a previous study of LANS— a [1^, 
we take a to be l/40th the box size. This choice was found to produce an optimal a— LES (in the sense of being 
optimal for the class of LANS— a models, with respect to the value of a). Finally, we review bounds on the size of 
the attractors and use these bounds to comment on the computational savings of the three regularizations viewed as 
SGS models. 



II. THE THREE REGULARIZATION MODELS 
A. Clark— a 

The incompressible Navier-Stokes equations are given in Cartesian coordinates by 

dtVi + dj{vjVi) + diV = v^jjv^, diVi = 0. (1) 

Filtering these equations with a convolution filter, L : z ^ z \n which z (resp. z) denotes the filtered (resp. unfiltered) 
field, yields 

dtUi + dj{ujUi) + diV + djfij = vdjjUi , (2) 

in which by convention we denote Ui = Vi and the Reynolds turbulence stress tensor, fij = ViVj — ViVj, represents 
the closure problem. Eq. ^ can represent either a LES or a SGS model. As the difference between the two is 
primarily philosophical (e.g., the scale at which filtering is applied, dissipative versus dispersive, the factor by which 
computational resolution may be decreased, etc.), we briefly deflne our terminology. Many LES include eddy-viscosity 
(i.e., djfij includes a VTdjjUi term such that vt ^ i'). This amounts to approximating the v = problem and no 
finite Reynolds number can be defined. More generally, a LES applies the filtering in the inertial range and reduces 
the necessary computational linear resolution by at least an order-of-magnitude. Different from this previous case, 
a SGS model employs a finite value of v (and a well-defined Reynolds number) and addresses instead the question: 
For a given Re, how far can we reduce the computational expense while retaining as much of the detailed large-scale 
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properties (such as the high-order statistics) as possible? For the case of LANS— a, for example, it has already been 
shown that the reduction in computational expense is rather modest (a factor of about 30). Therefore, while calling 
it a SGS model is justified by Eq. ([2]), the label LES does not really apply. 

It is the case for both LES and SGS models that though a single filtering is indicated in Eq. numerical solution 
implies a second filtering at the grid resolution (see, e.g., 0]). Systematic studies requirin g a database of computed 
solutions have been made in the past for LES [3§| and for the LANS— a regularization model [isl. [2ll. [45l| . These studies 
show that the ratio of the two filter widths (i.e., the sub-filter resolution) can affect greatly the model's performance. 
To avoid this complication, the sub-filter resolutions employed in this study are rather large. Determination of the 
optimal sub-filter resolution is a detailed study which should be undertaken for both Clark— a and Leray— a, but is 
beyond the scope of this present paper. 

It has been shown 0, [52 | that for all symmetric filters possessing a finite nonzero second moment, the first term of 
the reconstruction series for the turbulent sub-filter stress is 



dkUidkUj + ... (3) 

dk^ k=o 



where G{k) is the Fourier transform of the convolution kernel (G'(r) is the convolution kernel where [Lz]{r) = J G{r — 
r')z(r')d^r'). This approximation of the subgrid stress is then generic and is known as the Leonard tensor-diffusivity 
model [1^ (or, often, the Clark model [l^). Related to this model, the Clark— a model is (6| 

dfVi + Hdj{ujUi) -f diP + a^dj{dkUidkUj) = vdjjVi, (4) 

and its subgrid stress for a ^ 1 is given by 



T, 



n-^a^dkuAuj) = a^idkuAuj) + 0{a^). (5) 



Here the filter is the inverse of a Helmholtz operator, L = Ti ^ = (l — a^V^) ^. The Clark— a model conserves 
energy in the H\{u) norm instead of the Lp'{v) norm. 



dt 

in which the Clark— a energy Eq. is expressed as 



"^^^ = -2u^^ , (6) 



= ;^ ^ (u - V^u) • ud^a; = ^ j J^^ ■ ud^'x , (7) 
and the Clark— a energy dissipation rate is given by 

na = ^J^^uj- u,d^x , (8) 

where a; = V x v and a> = V x u. For L = H^^, we note that G = (1 + a'^k'^)^^ which implies that the turbulent 
sub-filter stress tensor for the tensor-diffusivity model given by Eq. ([3]) is 

f^^ =2a^{dkU^dkU,) + ... (9) 

which is proportional to the Clark— a stress tensor to second order in a. Hence, the a priori tests of [s^ should apply 
to Clark— a, at least at this order. 



1. Kdrmdn-Howarth equation for Clark— a 

In 1938, Karman and Howarth [l^ introduced the invariant theory of isotropic hydrodynamic turbulence, and 
derived from the Navier-Stokes equations the exact law relating the time derivative of the two-point velocity correlation 
to the divergence of the third-order correlation function. The corresponding Karman-Howarth theorem for LANS— a 
in the fluid case was derived in [2^. The relevance of the Karman-Howarth theorem for the study of turbulence 
cannot be underestimated. As a corollary, rigorous scaling laws in the inertial range can be deduced. In this section, 
we derive these results for the Clark— a case. 
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For the sake of simplicity, we consider the case i> = 0, since the dissipative terms may be added at any point in 
the derivation. We denote u' = u(x',i) and begin our investigation of the correlation dynamics by computing the 
ingredients of the partial derivative dt{viu'^). The Clark— a motion equation ([4|) may be rewritten as 

dtVi + dm{ViU,n + UiVm - UiUm + p5im " dnUidnUm) = 0. (10) 

Combining Eqs. ^ and ([5]), we arrive at the fluctuation- velocity equation, 

dt u'^ + <9„ (u; u\^+p'5j,^+a^G®Tf^)=0 , (11) 
where f^, = H^^a^r^ [L = H^^). Multiplying Eq. (jTU]) by u' and Eq. (fTT|) by Vi, then adding the result yields 



d d 

dt{Viu'A = {{VtU,n + UiVra - UiUra - Q;^9„Mj9„M„i)lt' ) + (pit'(5„„ - p'viSjrn) 

' OTm or., ■' 

d 



{{uW,n+0?G®Tf^^)v,), (12) 



9r„ 

where we have used statistical homogeneity 



d d d 

We symmetrize Eq. (|12p in the indices i, j by adding the corresponding equation for dt{vju'j). We then use homogeneity 
again as 

{vtu'jU^ + Vju'^u'„^) = - {v[ujUrn + v'jUiUm) , (14) 

and define the tensors 

Qf- = {v,u'^ + vju'i) , (15) 

(u.iU'j + Uju[) Vjn >, (16) 
= {{P'U-'j - P'Vj) Si„i + [pu', - p'vi) Sj„r) , (17) 
^fjm = < {dnUrdnU,n) u'j + (dnUjOnUm) u[ + G ® T^Jv, + 

G ® rfjv, > . (18) 

We can drop because the terms with the pressures p and p' vanish everywhere, as follows from the arguments 

of isotropy [IJ] . Finally, we obtain 



^* 2g = {T^^ ~ . (19) 



This is the Karman-Howarth equation for Clark— a (compare to Eq. (3.8) in [23 for LANS— a). 

By dimensional analysis, the energy dissipation rate in Clark— a is ~ dtQ^ ij and Eq. (jl9[) implies 

e^^-\{vu'+u^+'^u'). (20) 

For large scales {I ^ a), we recover the Navier-Stokes scaling known as the four-fifths law, < {Sv\\{l))^ >^ el 
Here, Sv^^ (l) = [v(x + 1) — v(x)] • l/l is the longitudinal increment of v. Strictly speaking, the four-fifths law expresses 
that the third-order longitudinal structure function of v, S'K/) = (((5u||)'^), is given in the inertial range in terms of 
the mean energy dissipation per unit mass e by 

^3 = -^el, (21) 
or, equivalently, that the fiux of energy across scales in the inertial range is constant. We also recover the Kolmogorov 
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1941 [3l|, m, [111 (hereafter, K41) energy spectrum, E{k)k ^ v'^ ^ e^/^P/'^ or, cquivalcntly, 

E{k) ~ e2/3fc-5/3. (22) 

For sub-filter scales (/ <S a), we have u ~ vl^ jo? and the first and third right-hand terms in Eq. (|20p are equivalent. 
In this case, we are left with two different possible scalings depending on the prefactors in Eq. (pO|) . If the first (or 
third) right-hand term is dominant, our scaling law becomes 

((<5«ll (0)^(^x^11 (0))-efZ- (23) 

For our sub-filter scale energy spectrum we would then have E^{k)k ~ uu ~ {e^)'^^^a'^^^ , or, cquivalcntly, 

E^ik) ^ (£^)2/3a2/3fc-i. (24) 



This result is the same as for the a— model |16[. If, however, the second right-hand term in Eq. ([20| is dominant, 
then the K41 results are recovered, with u substituted for v. In that case, one finds the alternative Clark— a sub-filter 
scale spectral energy scaling, 

E^ik)^k'/\ (25) 

2. Phenomenological arguments for Clark— a inertial range scaling 

We review here the derivation by dimensional analysis of the spectrum which follows the scaling ideas originally 
due to Kraichnan [3^ and which is developed more fully in Ref. [6|. In examining the nonlinear terms in Eq. ([4]), it 
is not entirely clear which of three possible scales for the average velocity for an eddy of size k~^, 

1/2 

= / \-^\'d'-) ^ (26) 



D 



1 



1/2 

Uk ■ Vferf^a; ) 

D 



or 



ui'^=iT^I \^^\'d'x] (28) 



1/2 

\uk\^d'x] 
Id J 

should result. Therefore, three corresponding "turnover times", t^, for such an eddy may be proposed 



4"^ - l/(fcC^i"') wfth n = 0, 1, 2 . (29) 

The term "turnover time" is used advisedly here, since only the velocity is composed of the fluid transport 
velocity. We define the (omnidirectional) spectral energy density, Ea{k)^ from the relation 

nOO n poo 

Ea= I (f)Ea{k)dadk^ / Eo,{k)dk . 
Jo J Jo 

Since, ■ = • Vfe/(1 + a^fc^) = Ea{k)/{1 + a^k'^), we have 

1/2 



(30) 



( / £;a(fc)(l + a2fc2)(i-")dfcj ~ (fc£:„(fc)(l + a2fc2)(i-«)^ (31) 
Then, the total energy dissipation rate, ej^ is related to the spectral energy density by 

So - (4"')"' / E^{k)dk ~ fc2t/(")i^„(fc) ^ fc5/2^C(fc)3/2(i ^ a2fc2)(l-„)/2^ (33) 
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which yields, finaUy, the predicted energy spectra for Clark— a, E^{k), 

E^ik) - (e^i )'/3fc-^/3(l + a2fc2)("-i)/3. (33) 
For scales much larger than a {ak <^ 1) the Kolmogorov scaling for Navier-stokes is recovered, 

E^{k) ^ {e^^f"k-''\ (34) 
whereas for scales much smaller {ak ^ 1), the spectrum becomes 

E^ik) ~ (eC)2/3^(2(„-l)/3^(2„-7)/3_ (35) 

These arguments constrain the Clark— a sub-filter scale spectrum to lie between k~^ and k~'^^^. 

B. Leray— a 

The Leray model in Cartesian coordinates is 

dfVi + dj{ujVi) + diP = vdjjUi diVi = 0, (36) 



where the flow is advected by a smoothed velocity, u. By comparison with Eq. ^ we see that the Leray model 



approximates the subgrid stress as r,^ ~ L{ujVi) — UjUi, or, with L = Ti ^, 



H a [dkUidkUj + dkUidjUk). (37) 



As has been noted previously [2^, the subgrid stress of Clark— a in Eq. ([51) is a truncation of the subgrid stress of 
Leray— a, in Eq. (|37p . For Leray— a, the L'^{v) norm is the quadratic invariant that is identified with energy, 

dE 

^ = -2un, (38) 



where 



and 



i=4/J|v|Vx. ,39) 



(1 = ^ J^^l^l'd'^. (40) 

As was pointed out in [sij, the incompressibility of the velocity field v only implies a divergenceless filtered velocity 
u under certain boundary conditions for Leray— a. When diUi ^ 0, the energy E = ^ Jjj |vp is no longer conserved 
(helicity and Kelvin's theorem are not conserved for Leray— a). In our numerical study, we employ periodic boundary 
conditions, for which diVi = implies diUi = and Leray— a conserves energy in the usual sense of L^{v). 

1. Kdrmdn-Howarth equation for Leray— a 

In this section we derive the Karman-Howarth equation for the Leray— a case. Following Section fll A 11 we begin 
our investigation of the correlation dynamics by computing the ingredients of the partial derivative dt{viv'^). Eq. ([36]) 
may be rewritten as: 

dtVi+ dm{v^Um + P5tm) (41) 

Multiplying Eq. ([JT]) by v'j yields 

dMv^) = J—{v^u„y^) + —{Pv'^5,m)- (42) 
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Wc can make this equation symmetric in the indices i, j, by adding the equation for dt{vjv'^). We define the tensors 

Qf^ = {v,v'^+v,v';), (43) 
%jm = {{v^v'j + Vjv';) u^,) , (44) 

nf,™ = {Pv'^S„n + PvlS.m) . (45) 

Again, we may drop n^„j because the terms with the pressure P vanish everywhere and thereby obtain 

= 5^^'™- (46) 

This is the Karman-Howarth equation for Leray— a. 

The energy dissipation rate for Leray— a is denoted by and it satisfies ~ dtQ^ij. By dimensional analysis, 
Eq. (j46|) implies 

£^ ~ jv^u. (47) 

For large scales (l ^ a), we recover the Navier-Stokes scaling, Eqs. (|2T|) and (|22|) . For sub-filter scales (/ <^ a) our 
scaling law becomes 

((Jz;||(0)^^«||(/)))~e^Z. (48) 

For our small scale energy spectrum we would then have E^{k)k ^ v"^ ~ (£^)^/^a"'/^fc^/'^ (where we employed 
u ~ vP/a'^), or, equivalently, cf. Eq. (|25|) . 

i5^(fc)-(e^)2/3a4/3^-i/3^ (49) 

2. Phenomenological arguments for Leray— a mertial range scaling 

We review here the derivation by dimensional analysis of the spectrum for Leray— a as we did for Clark— a in 
Section III A 21 This analysis is developed more fully in Ref. [l^l ■ We argue again that there are three possible scales 

for the average velocity for an eddy of size k~^, Eqs. (pS)) . (P7|) . and (pS)) . with the turn-over time, t^"^ given by Eq. 
(UHl). Since, ul ^vl/{l + a^k^f ^ E{k) / [l + a^k'^f , we have 

/ r \ 



Uk ~ + a k^ydkj ~ (fc£:(fc)(l -f a^k^y) ' (50) 

Then, the total energy dissipation rate, is related to the spectral energy density by 

£L ^ (t(^»))-l f E{k)dk - eu\!'^E{k) ~ fc5/2i;(fc)3/2(l + a2fc2)-"/2^ (5I) 



which yields, finally, the predicted energy spectra for Leray— a, E^{k), 

E^{k) ~ (e-^)2/3fc-5/3(i + a2fc2)„/3_ (52) 

For scales much larger than a {ak ^ 1) the K41 spectrum is recovered, Eq. (P^ . and for scales much smaller {ak 1) 
the spectrum is 

E'^ik) ~ (£i)2/3^2„/3^(2„-5)/3_ (53) 

These arguments constrain the Leray— a sub-filter scale spectrum to lie between k~^^^ and k~^^^. 
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C. LANS-a 



LANS— a is given by 

dtVi + dj{ujVi) + diir + vjdiUj ~ vdjjVi, diVi = 0, (54) 
For LANS— a, the usual choice of filter is again L = Tir^ With this filter, the subgrid stress tensor is given by 

fj" = H^^a'^idrnUidmUj + drnUidjUm ~ diUmdjUm.)- (55) 

As has been previously noted [1^, the subgrid stress of Leray— a, Eq. ([37)) . is a truncation of the subgrid stress of 
LANS— a Eq. ((55|) . Like Clark— a, energy is conserved in the H^{u) norm instead of the L'^{v) norm. Additionally, 
LANS— a is the only model of the three examined here that conserves a form of the helicity (and Kelvin's circulation 
theorem). 

For LANS— a in the fluid case the Karman-Howarth theorem was derived in [23 | . We summarize here the dimensional 
analysis argument for the LANS— a inertial range scaling that follows from this theorem, beginning from Equation 
(3.8) in Ref. [1^. In the statistically isotropic and homogeneous case, without external forces and with v — 0, taking 
the dot product of Eq. ((54|) with u'j yields the equation 



tjm/ 



(56) 



The trace of this equation is the Fourier transform of the detailed energy balance for LANS— a; 

Q?j^{v,u^+v,u,) (57) 

is the second-order correlation tensor while 

"^jm = {{ViU^ + VjUi + V^Uj + VjUi)Um) , (58) 



and 

S"jm = (^{drnUidiUi)uj + {d„iUidjUi)ui + (G (g) T"J)vi + (G (g) r^„')wj) , (59) 

are the third-order correlation tensors for LANS— a and f°' = H^^a^rf' is the sub-filter scale stress tensor. For a ~ 

''J ''J 

this reduces to the well-known relation derived by Karman and Howarth. The energy dissipation rate for LANS— a, 
Eq,, satisfies £„ oc dtQ^ij- By dimensional analysis in Eq. (ISS)) we arrive at 

1 

-(vu^ + —u^). (60) 

For large scales (/ 3> a), we recover the Navier-Stokes scaling Eqs. (PT|) and (P^ . For sub-filter scales (/ ^ a) our 
scaling law becomes Eq. (j23p and our sub-filter scale spectra is given by 

E^{k) ^ el/'a'/'k-\ (61) 

In this case, by the phenomenological arguments, we know that eddies of size are advected by the smoothed 
velocity, Eq. ((58)) . This scaling is confirmed in Ref. [i^ but it coexists with a energy spectrum corresponding to 
"enslaved rigid bodies" or "polymerized" portions of fiuid which do not contribute to the turbulent energy cascade. 



III. SUB-FILTER-SCALE PHYSICS 



Only by examining the sub-filter scales can we hope to derive new, improved models, and, ultimately, to gain an 
understanding of turbulence. A knowledge of the differences between closures and Navier-Stokes is fundamental to 
enable the derivation of better physical models of turbulence at small scales. In this section, then, we will be interested 
in both the similarities and the differences between the regularizations and Navier-Stokes. A more immediate goal of 
predicting the computational savings at higher Reynolds numbers can be achieved through the correct prediction of 
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the scahng at small scales. 

To this end, we compute numerical solutions to Eqs. ((1]), (15^. a nd ([51)1 in a three-dimensional (3D) cube with 
periodic boundary conditions using a parallel pseudospectral code [22l. [23l | . We employ a Taylor-Green forcing , 



F 



sin k()X cos kgy cos koz 
— cos fcox sin fcoy cos fco^^ 




(62) 



(with kQ ~ 2), and employ dynamic control |41| to maintain a nearly constant energy with time. The Taylor-Green 
forcing, Eq. ([62|) . is not a solution of the Eulcr's equations, and as a result small scales are generated rapidly. The 
resulting flow models the fluid between counter-rotating cylinders @ and it has been widely used to study turbulence, 
including studies in the context of the generation of magnetic flelds through dynamo instability [i^. We define 

1/2 

the Taylor microscale as A ~ 2ii ^J{v^)/iJJ^ ^ and the mean velocity fluctuation as Vr^ns = {2j^E{k)dky'\ The 
Taylor microscale Reynolds number is defined by R\ = Vrms-^/^ s-nd the Reynolds number based on a unit length is 



Re 



The Clark— a, Leray— a, and LANS— a equations (as well as other SGS models based on spectral filters) are easy 
to implement in spectral or pseudospectral methods. As an example, in Fourier based pseudospectral methods, the 
Helmholtz differential operator can be inverted to obtain 7i~^(fc) = (1 + a^fc^)~^, where the hat denotes Fourier 
transformed. In this way, the filter reduces to an algebraic operation and Eqs. Q, and (|54p can be solved 

numerically at almost no extra cost. If other numerical methods are used, the inversion can be circumvented for 
example by expanding the inverse of the Helmholtz operator into higher orders of the Laplacian operator [13, [S^] • 

To compare the three regularizations (Clark— a, Leray— a, and LANS— a) we compute a fully resolved DNS of the 
Navier-Stokes equations at a resolution of 1024^ (y = Z ■ 10^^, Re « 3300) and model runs with the exact same 
conditions at a resolution of 384^. The details of the flow dynamics of the DNS have already been given [l], |40| . 
In particular, the Reynolds number based on the integral scale £ = 27r / E{k)k^^dk/E w 1.2 (where E is the total 
energy) is Res, = UZ/v Ki 3900, where U is the r.m.s. velocity and the Reynolds number based on the Taylor scale is 
Rx ~ 790. The DNS was run for nine turnover times {Z/U) (in the following results, time t is in units of the turnover 
time). We employ a filter width of a = 27r/13 for which LANS— a exhibits both Navier-Stokes and LANS— a inertial 
ranges in the third-order structure function [45l |. From these we hope to obtain the behavior of the models for scales 
much smaller than a 




FIG. 1: (Color online) Time evolution of the enstrophy {{i^^) for Leray— a and in the DNS, (tu • ZJ) for LANS— a and Clark— a). 
DNS {Re « 3300) is shown as solid black lines, LANS— a as dotted red, Clark— a as dashed green, and Leray— a as blue dash- 
dotted. The cyan long-dashed line represents a 512"^ DNS {Re « 1300, R\ ~ 490). Here each run is calculated only until it 
reaches a statistical steady state. Leray- a reduces the dissipation, e = u {^^), and increases the time scale to reach a statistical 
turbulent steady-state. Both effects are greater as a is increased. By comparison with the Re ~ 1300 run, we see that these 
two effects are consistent with a reduced effective Reynolds number. A smaller reduction in flux (but not an increase in time 
to steady-state) is also observed for LANS— a and is likely related to its rigid bodies. 

In Fig. [l]we present the time evolution of the enstrophy ((w^) for Leray— a and DNS, (u) - TD) for LANS— a and 
Clark— a), which is proportional to the dissipation {e = v (w^) or = i/ (u) • a;) depending on the case). Also shown 
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is a wcU-rcsolvcd 512'^ DNS of a less turbulent flow {ly ~ 1.5 • lO"'^, Re « 1300, R\ w 490, (cyan online) long-dashed 
line). Here each run is calculated only until it reaches a statistically steady state. We see that the dissipation is greatly 
reduced and the time-scale to reach a statistically steady state is increased for Leray— a. We see by comparison with 
the Re « 1300 DNS that this reduced dissipation could result from a reduced effective Reynolds number in the 
Leray— a run. For LANS— a, the dissipation is decreased although the time to reach steady state is not increased. 
This is probably related to the enslavement of its rigid body regions which would have no internal dissipation. Of 
the three models, Clark— a most resembles the total dissipation for a large range of a. Indeed, as it is the order 
approximation of Navier-Stokes, this dissipation behavior for Clark— a may continue to hold until a becomes quite 
large. 




FIG. 2: (Color online) Spectra compensated by K41 for 1024^ DNS {Re f» 3300) averaged over t = [8.25,9]. Labels are as in 
Fig. [1] The vertical dashed line indicates ka = 27r/a. Compensated spectra for 384"^ LANS— a averaged over t G [10.8, 11.6], for 
384=* Clark-a over t € [11.8, 12.6], for 384=* Leray-a over t G [26.3, 27], and for the Re ^ 1300 DNS over t G [18.1, 18.9]. Due 
to the large disparity in times to reach a turbulent steady-state, the time intervals chosen to average over also differ greatly. 
Clark— a approximates the predicted spectrum, Eq. (|24|l . and not k^^^, Eq. (|25|l . nor another possible spectrum, Eq. (|33p . 
The spectrum of Leray-a is very similar (for k G [5, 20]) to that of the Re ~ 1300 DNS. The positive but shallower than k^ 
LANS— a spectrum observed here has previously been reported (see text). The spectra of all three regularizations clearly differ 
from that of Navier-Stokes. 



Figure [2] shows a comparison of the energy spectrum at the turbulent steady state, for all the runs in Fig. [T] The 
isotropic energy spectra are calculated as follows, 



fccff<fc+T 

m = 

Ka>k-i 



W^(fccff) + Wy(fccff) + vlikcfs) 



(63) 



where k^s — + ky + (the Hf{u) norm is employed for Clark— a and LANS— a). The length scale a is indicated 

by a vertical dashed line and the plotted energy spectra are compensated by k^^^ (i.e., leading to a flat K41 
spectrum). The energy flux in the DNS is constant in a wide range of scales, but the compensated spectrum has a 
more complex structure. The salient features of this spectrum are well-known from previous studies [30|. Small scales 
before the dissipative range show the so-called bottleneck effect with a slope shallower than k~^/^. On the other hand, 
larger scales have a tendency to develop a spectrum slightly steeper than k~^/'^ because of intermittency corrections, 
an effect that becomes clear in the simulation performed at larger spatial resolution on a grid of 4096^ points [sO ]. 
From Fig. [5] it is clear that the Clark— a spectral behavior is close to the predicted fc~^ spectrum, rather than the 
^1/3 from Eq. (|25p . or the other possible spectrum from Eq. ([33]) . Likewise, the positive LANS— a spectrum 
observed here approaches k^ with increasing resolution as the sub-filter scales are fully resolved [1^. Leray-a, on 
the other hand, possesses a very steep sub-filter scale spectrum as well as enhanced large-scale energy as has been 
previously observed [2l| . The results indicate that solutions to Leray-a are the most strongly regularized of the three 
regularizations. 

The spectrum of Leray-a in Fig. [2] gives a good approximation to the Re w 1300 DNS in the range k G [5, 20] (i.e., 
io V = 1.5 • 10~'^ rather than to = 3 • 10~^ which was employed). This result, Leray— a's increased characteristic 
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time scales and its reduced dissipation imply that the Leray— a model is operating at a much lower effective Reynolds 
number. This is also clear from the rapid drop in the spectrum at small scales, shown in Fig. [51 Indeed, we can build 
an effective Reynolds number in the large scales as R^s — ^^'^H'^l'^ jv. Since £ is controlled in this simulation by the 
forcing scale, the drop in the dissipation rate implies a reduced nonlinearity in Leray— a. This is also consistent with a 
direct comparison of the nonlinear terms in Leray— a with, for instance, LANS— a. The nonlinear terms in LANS— a, 
Eq. ([54]) . may be written as u ■ Vv + Vu-'" • v (where the sufRx T denotes a transposition), while the nonlinear term 
in Leray— a, Eq. is only u • Vv. Both nonlinear terms in LANS— a are of order 0(1); so the absence of one of 

the nonlinear terms in Leray— a could be understood as a reduction in the nonlinearity. 




FIG. 3: (Color online) Third-order structure function associated with the Karman-Howarth equation, L'"'^'(Z), versus length, 
I. Labels are as in Fig. [T] The vertical dashed lines indicate the length a. The Clark— a result is consistent with a m^u ~ 1} 
scaling, Eq. (|23|l . and clearly inconsistent with a m^v ~ (y? I) scaling as would arise from the middle term in Eq. (I20|l . 
The results for Leray— q are again consistent with a reduced effective Re. 

Validation of the Karman-Howarth equation scalings, Eqs. pS)) and ((48)) . enables us to measure scaling laws in the 
inertial range and, thus, compare the intermittency properties of the models. The third-order correlations involved 
in the theorems, namely 

L"(0^((<5^,l(Of('5t'll(0)) (64) 

for Clark— a and LANS— a, and 

L^0^(('5^;||(0)'('^^||(0)) (65) 

for Leray— a, and L{1) = S'!^{1) for Navier-Stokes, are plotted versus / in Fig. [3l In Fig. [3] we can see validation 
of the Karman-Howarth scaling for scales smaller than a for both LANS— a and Clark— a. In particular, we note 
the observed scaling for Clark— a verifies the vu^ ~ I scaling and not the (theoretically possible) vu^ ~ {u^ ~ I) 
scaling. The predicted scaling is not observed in Leray due to its reduced effective Reynolds number. With these 
scalings in hand, we may proceed to observe the scaling of the longitudinal structure functions, 

S;il)^{{Sv\\'r/'}, (66) 

where we again replace the norm, | |), for the norm, {{5v^^)^), in the case of Clark— a and LANS— a. 
We utilize the extended self-similarity (ESS) hypothesis [1, d, which proposes the scaling 

5;(/)cx(l("^^)(0)^'' , 

and normalize the results by ^3 to better visualize the deviation from linearity (which serves as a measure of intermit- 
tency). As we will show in the next section, our flow is anisotropic in the z— direction. Therefore, structure hmctions 
are computed in horizontal planes only. The results are displayed in Fig. (¥] 

In Fig. m we may observe the intermittency properties of the models at sub-filter scales. We note a reduced 
intermittency for both Leray— a and the Re w 1300 DNS. This is consistent with the smoother, more laminar fields 
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FIG. 4: (Color online) Normalized structure function scaling exponent £,p/(,3 versus order p. The dashed line indicates K41 
scaling and the solid line the She-Leveque formula |4^. The DNS results are indicated by black X's, LANS— « by red asterisks, 
Clark— a by green diamonds, Leray— a by blue triangles, and the Re ~ 1300 DNS results are shown by cyan pluses. Leray— q 
is less intermittent consistent with the smoother field produced by a lower Re flow. Clark— a is more intermittent than Navier- 
Stokes at sub-filter scales. LANS— a is less intermittent than Clark— a, likely due to the influence of its rigid bodies (see 
text). 

(due to the reduction of the eflFective Re) possessed by both. Interestingly, though LANS— a and Clark— a both possess 
the same cascade scaling (Eq. as confirmed in Fig. [3]), the Clark— a model is markedly more intermittent than 

LANS— a. If artificially truncated local interactions (in spectral space) is taken as a cause of enhanced intermittency 
[isl . [ssj . then the increased intermittency observed in Clark— a is the expected result of truncation of the higher- 
order terms in the sub-filter stress tensor. Moreover, if the LANS— a's ~ spectrum is indeed associated with rigid 
bodies, these would serve to decrease the intermittency (no internal degrees of freedom being available in a rigid body) 
which is consistent with the results shown here. Due to this effect, LANS— a of the three regularization models most 
resembles the high-order intermittency of Navier-Stokes at sub-filter scales. 

IV. SGS POTENTIAL OF THE REGULARIZATIONS 
A. Reproduction of super-filter scale properties 

The differences at sub-filter scales between the rcgularizations and Navier-Stokes are important to understand how 
the models may be improved upon. From a practical standpoint, an equally important question is how they predict 
the supcr-filter-scale properties of a DNS when employed as models. This gives an indication of their SGS modeling 
potential. For this we choose a = 27r/40 corresponding to an optimal a— LES [i^. Note that the value of a has been 
optimized for neither Clark— a nor Leray— a; as a consequence, these models might perform better in other parameter 
regimes than the results indicate in this study. 

Fig. Ogives the time evolution of the enstrophy of the DNS and the models along with that of an under-resolved 
Navier-Stokes solution at a resolution of 384^ (i^ = 3 • 10^"*^, (pink online) dash-triple-dotted line). We see that both 
LANS— a and Clark— a reproduce the proper amount of dissipation and are within 10% of the time required by the 
DNS to reach a statistical turbulent steady state. As has been observed before, Leray— a is under-dissipative [2l[. 
We also note that it takes longer than the other models to reach a steady state even with the smaller filter width 
(27r/40 as opposed to 27r/13). When compared to the larger a case, we see that the dissipation is much greater and 
the time-scale to reach a turbulent steady state is decreased for Leray— a. 

Compensated spectra averaged over several eddy turn-over times are shown for the SGS case (i.e., ka = 40) in Fig. 
|6l Note that as the subgrid models are averaged over a different time interval, no meaningful comparison to the DNS 
is possible for k < kp ^ 3. Even without an optimal choice for the value of a, Clark— a best reproduces the DNS 
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FIG. 5: (Color online) Time evolution of the enstrophy ((<^^) for Leray— a and in the DNS, {ui ■ uJ) for LANS— a and Clark— q). 
DNS {Re ~ 3300) is shown as solid black lines, LANS— a as dotted red, Clark— a as dashed green, and Leray— a as blue dash- 
dotted. An under-resolved (384^) Navier-Stokes run is shown as a pink dash-triple-dotted line. 




FIG. 6: (Color online) Spectra compensated by K41 for 1024^ DNS {Re « 3300) averaged over t — [8.25, 9]. Labels are as in 
Fig. [5l The wavenumber corresponding to the filter width, ka , is shown as a vertical dashed line. 384^ simulations are averaged 
over t £ [15,20]. Note that to make a comparison for most wavenumbers, the spectra must be averaged within a turbulent 
steady-state. Therefore, as the subgrid models are averaged over a different time interval, there is no meaningful comparison 
to the DNS for < 3. For LANS— a we observe a contamination of the super-filter-scale spectrum (at k G [9,30]) related 
to the steep sub-filter-scale spectrum and the conservation of energy. Even though a different a (optimized with respect to 
spectral energy prediction at this numerical resolution) may provide better results for Clark— a, this model does very well at 
reproducing the large-scale energy spectrum. Leray— q's performance is the poorest. 

spectrum for scales larger than a. Wc compute root-mcan-square spectral errors as recently introduced in Rcf. [3^ : 

Etu. pnE„.o,edk) - E{k)r y^' 

where kp is the wavenumber for the forcing scale, E{k) is the DNS spectrum (in the LF'{v) norm), and Emodei{k) is 
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the subgrid model spectrum (in the appropriate norm). Another measure introduced in Ref. [39| is given by 



With p = 0, we find the error in the total energy, = ee- As this is dynamically controlled in our experiment, we 
find zero in all cases. For p = 2, we find the error in the total dissipation, £3 = eg, which is observed in Fig. [5l 
Every deviation from the DNS spectrum is counted positive, however, in e^. For p = we find the error in the energy 
spectrum: in decreasing order, e(5 = 0.24 for Leray-a, 0.23 for the under-resolved 3843, 9.20 for LANS-a, and 0.16 
for Clark— a. Both LANS— a and Clark— a improve the estimate over the under-resolved run, but Clark— a makes the 
best prediction. We see that only Clark— a improves the estimate of the power spectrum at this resolution for each 
scale considered separately (see Fig. Leray— a performs the poorest of the three regularization models, but it is 
also not optimized. As previously argued, its effective Reynolds number is too low to accurately model the DNS flow. 
Either a decrease in the viscosity z/, or a decrease in the filter size a (and, hence, an increase in the nonlinearity) , or 
both would likely improve the accuracy of Leray-a as an SGS model. Due to its frozen-in (or enslaved) rigid-body 
regions and its conservation of total energy, the LANS— a model cannot reproduce the DNS spectrum at super-filter 
scales unless a is only a few times larger than the dissipation scale [45l |. 

Another measure of the success of a subgrid model is the reproduction of structures in the fiow. In Figure[7]we have 
3D volume rendering of the enstrophy density {lj-uj for LANS— a and Clark— a) for the DNS, the three SGS-model 
simulations {ka = 40), the 384"^ under-resolved Navicr-Stokcs solution, all at a Reynolds number of « 3300, and the 
Re w 1300 DNS. Due to the late times depicted (longer than a Lyapunov time) there can be no point-by-point 
comparison between the simulations. Instead, we note that there are four horizontal bands where the forcing causes 
a maximum shear. This large-scale feature of the flow is missing only from Leray-a and the Re ~ 1300 run. The 
three other runs reproduce this feature well (note that the apparently thicker tubes present in Clark— a are vortex 
tube mergers). The results lead again to the conclusion that the under-resolved Navier-Stokes, the Clark— a, and the 
LANS— a models are better subgrid models than Leray-a due to its reduced effective Re. 

For the SGS models, the predicted from the Karman-Howarth theorem for Navicr-Stokcs is well-reproduced by 
all models at super-filter scales (not depicted here). We may then proceed in Fig. [8] to analyze the SGS model 
intcrmittency results. We see that all models reproduce the intermittency up to the tenth-order moment within the 
error bars (although there is a small decrease in intermittency for Leray-a). Thus, we conclude that with adequately 
chosen values of a (and of i> for Leray-a), all three models can reproduce the intermittency of the DNS (to within 
the error bars). 

The sub-filter-scale physics of Leray-a shows that it possesses the smoothest solutions of the three models and 
reduces the effective Re. We have seen that this strongly hampers its effectiveness as a SGS model. "Rigid bodies" are 
observed in the sub-filter scales of LANS— a that strongly influence even the super-filter-scale energy spectrum 
but not the sub-filter-scale dissipation nor intermittency properties. These affects also carry over to its application 
as a SGS model in that very small filter widths are required to properly predict the large-scale spectrum. Clark— a's 
approximately sub-filter energy spectrum is the closest to k~^^^ of the three models and is seen to cause the least 
contamination of the super-filter-scale spectrum when employed as a SGS model. Finally, when the filter width is 
small enough, the enhanced intermittency of Clark— a is nearly eliminated. 



The rationale behind using a SGS model is that it leads to adequate solutions at a reduced computational cost, 
because it computes fewer dof; indeed, for an SGS model, the ratio of Navier-Stokes's dof to the model's dof, a 
prediction for memory savings and hence computation time savings for numerical simulation, is a crucial factor. 
Consequently, analytical bounds on the sizes of the attractors for the three regularization subgrid models may be 
useful indicators of their computational savings. The dof for LANS— a is derived in [l^ and confirmed in (45^ . 




(68) 



B. Computational gains 



dof^ oc -i?e3/^ 



(69) 



(a) Re « 3300 DNS 



(b) Clark-a 





(c) Leray— a 



(d) LANS-a 





(e) Under-resolved Navier-Stokes (f) Re « 1300 DNS 

FIG. 7: (Color online) Volume rendering of the enstrophy density uj^ {u> ■ u) for LANS— a and Clark— q). The four lengths 
depicted are integral length scale £, Taylor scale A, filter width a, and dissipative scale 7^k as calculated separately for each 
simulation. For (a) Re f» 3300 DNS, (b) Clark— a, (d) LANS— a, and (e) under-resolved Navier-Stokes the snapshot is for 
t — 9. For (c) Leray— a it is for t = 16 and for (d) Re ~ 1300 DNS it for t = 19 corresponding to their slower development of 
turbulence. For Leray— a the location of vortex tubes are consistent with a lower Re flow while the other models (including 
under-resolving) reproduce the large-scale pattern of the flow well. The color scale indicates the strength of the enstrophy 
density with purple stronger than green. 



where L is the integral scale (or domain size). We may compare this to the dof for Navier Stokes, 

3 



dof. 



NS 



L 



Re^/\ 



(70) 



which immediately yields 



dof 



NS 



dof a 



q)Re^'\ 



(71) 



It was found, however, that to reproduce the super-filter-scale energy spectrum of an equivalent DNS, the filter-width 



a must be no larger than a few times the dissipation scale, "qx |45j . This is the result of the "polymerization" of 
the flow in LANS— a, and the associated E{k) ^ scaling at sub-filter scales and consequent contamination at 
super-filter scales via energy conservation. With this added caveat, it follows that the reduction in dof is independent 
of Re (and a net factor of about 10). Our study here illustrates that the high-order structure functions may be 
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FIG. 8: (Color online) Normalized structure function scaling exponent £,p/(,3 versus order p. The dashed line indicates K41 
scaling and the solid line the She-Leveque formula [4^. The DNS results are indicated by black X's. LANS— « by red asterisks, 
Clark— a by green diamonds, Leray— a by blue triangles, and pink boxes for the under-resolved Navier-Stokes run. With a 
small enough filter-width, a, the intermittency properties of the DNS can be reproduced with all three models. 



reproduced for much larger values of a. Therefore, in applications where the spectrum is not of great concern, much 
greater reduction in numerical resolution would be feasible. 

For Clark— a there is an upper bound on the Hausdorff, dn, and fractal, dp, dimensions of the attractor. 



a 



dH<dF<C{^] [-] , (72) 



where 77^ is the Kolmogorov dissipation length scale corresponding to the Clark— a model @. From its observed 
spectrum, we may estimate r]j^ or, equivalently, ^ Xjrf^. For dissipation the large wavenumbers dominate and, 
therefore, combining the Clark— a energy balance Eq. ^ with its sub-filter scale energy spectrum Eq. ([24]) allows us 
to implicitly specify its dissipation wavcnumber, /c^, by 

eE'^{k)dk ^ P(e^)2/3«2/3fc^lrffc ^ (^.C)2/3^2/3(^C)2^ (73) 



Then we have. 



It follows that 



S j,l/2„l/3- ^'^1 



dof^s ^ j^e^/4(")3/4^-i_ (75) 

do J Clark L 

This is similar to the prediction for LANS— a, but as energy spectra are more easily reproduced for larger values of a 
than with LANS— a (but not the intermittency properties), it may be the case that a is not tied to the Kolmogorov 
dissipation scale rjK- If so, then the computational saving might increase as Re^^^ which is promising for use of 
Clark— a as an LES model. This conclusion is bolstered to the extent that the results in Section IIV Al for ka = 40 
(a « 7r]K) arc acceptable. If even further separation from the dissipative scale is not possible, there is still a greater 
reduction in dof (a factor of 20) for Clark— a than for LANS— a. 

For Leray— a, we have the following upper boimds on the Hausdorff dimension, dn, and fractal dimension, dp, of 
the global attractor, 



d,<dp<c[-^) + (76) 
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where is the dissipation length scale for Leray— a p^ . Again, we estimate the dissipation wavenumber for Leray— a 
kfi - 1/Vk- Then, fromEqs. ^ and (gHl), that is, assuming the k spectrum resulting from the Karman-Howarth 
equation, we find 



Consequently we have. 



It follows that 



k^E'^{k)dk ~ / " fc2(£i)2/3^4/3^-l/3^^ _ (£^)2/3a4/3(fci)8/3. (77) 



3/8 I/,- C^^) 



(79) 



Our results suggest that for an effective LES the viscosity vl must be chosen to be smaller than ly. This leads to an 
upper bound on the computational savings for Leray— a, 



dof Leray (1 + ^f'^^ 

If we further assume that a is directly proportional to the dissipative scale rjK, we arrive at 



'^"^^^ <Ci?e27/56 (81) 



dof 

Leray 

which is not exceedingly promising for use as a LES. All such estimates are, however, purely conjectural until the 
proper choices of a and vl are determined. 



V. DISCUSSION 



We derived the Karman-Howarth equations for the Leray— a and Clark— a models. These two models may be 
viewed as successive truncations of the sub-filter scale stress of the Lagrangian-Averaged Navicr-Stokcs a— model 
(LANS— a). In the case of Clark— a two different inertial range scalings follow from the dimensional analysis of this 
equation. The case of Leray— a is simpler as a single scaling is predicted. This is the case for Navier-Stokes and 
LANS— a as well. To our knowledge, we computed the first numerical solution of the Clark— a model, the results of 
which are encouraging for further study. We compared these to solutions for a 1024^ DNS under periodic boundary 
conditions {v = i - lO""*, Re k, 3300) using a 384^ resolution under the same exact conditions for LANS— a, Leray— a, 
Clark— a, and an under-resolved 384^ solution of the Navier-Stokes equations. We employed two different filter widths 
a. The first choice a = 27r/13 was used to understand the sub-filter-scale physics and the second choice a = 27r/40 
was employed to test the SGS potential of the models. In comparing these two choices, we found for Leray— a that 
an increase in a substantially decreases the nonlinearity (and hence decreases the effective Reynolds number Re). 
For this reason, we were unable to confirm either the inertial range scaling from its Karman-Howarth equation or 
its sub-filter scale energy spectrum. For Clark— a we were able to determine the dominant Karman-Howarth inertial 
range scaling to be v?v ^ I which leads to the associated energy spectrum, also indicated by our results. 

The performance of the three regularizations as SGS models (for a resolution of 384^ and = 40) was comparable 
to that of the undcr-resolvcd Navicr-Stokcs solution in reproducing the DNS energy spectrum at super-filter scales. 
Only Clark— a showed a clear improvement in approximating the spectrum. From 3D volume rendering of enstrophy 
density we found that Clark— a and LANS— a were comparable to the under-resolved solution. Even at a = 27r/40, 
Leray— a's 3D spatial structures are consistent with a significantly reduced Re flow (e.g., comparable to a Re w 1300 
DNS). We note that the value of a was chosen optimally for LANS— a at the resolution of 384^, and that for Clark— a 
(and especially for Leray— a) smaller resolutions (greater computational savings) may have comparable results for this 
value of a. Such a comparison is beyond the scope of the present work. 

Although LANS— a and Clark— a exhibit the same inertial range scaling arising from similarities in their Karman- 
Howarth equations, Clark— a is decidedly more intermittent than Navier-Stokes at sub-filter scales. At the same 
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time, LANS— a is only slightly more intermittent than Navicr-Stokes. These results are consistent with the artificial 
tnmcation of local nonlinear interactions (in spectral space) in the SGS stress tensor of each model. This effect is 
reduced for LANS— a by the "rigid-body regions" enslaved in its larger scale flow which possess no internal degrees of 
freedom. The reduced intermittency observed for Leray— a is related to its smoother, more laminar fields as a result 
of its reduced effective Re. 

Finally, we analyzed the reduction in the number of dof in the models, as compared to Navier-Stokes (and, hence, 
their LES potential based on their computational savings). We noted that as LANS— a reproduces the intermittency 
properties of a DNS quite well even for larger values of a, some further reduction in numerical saving might be achieved 
provided the contamination due to its rigid-body energy spectrum were not important in a given application. As 
Clark— a possesses a similar reduction in dof to LANS— a, its LES potential is tied to the optimal value of a for 
LES. Our study indicates that Clark— a may be applicable (especially with regards to the energy spectrum) for larger 
values of a than LANS— a. In fact, if its optimal value is not a function of Re, the computational resolution savings 
increases as Re^/^ for Clark-a. For the case of Leray— a, the prediction is complicated by the effective reduction in 
Re as a increases. Prediction of optimized values of a and of effective dissipation i^l are required to assess its LES 
potential. Future work should include such a study for both Leray— a and Clark— a. 

All three regularizations were shown to be successful, in that their control of the flow gradient reduces the degrees 
of freedom and saves computation while preserving a properly defined Reynolds number (albeit for Leray— a that 
definition is not yet demonstrated). Clark— a accurately reproduces the total dissipation, the time scale to obtain a 
turbulent statistical steady-state and the large-scale energy spectrum of a DNS. These results seem to result from 
Clark— a being an order approximation of Navier-Stokes. We have shown that Leray— a reduces the effective 
Reynolds number of the flow. The last of the three models, LANS— a restores Kelvin's circulation theorem (advected 
by a smoothed velocity) and the conservation of a form of hclicity. Using spectra as a measure of the success of a 
subgrid model, LANS— a is less than optimal, due to its contamination of the super-filter-scale spectrum. However, 
other measures of the success of a subgrid model are possible: for example, in regard to intermittency. LANS— a may 
be considered a superior model. For Clark— a, intermittency may be a function of filter width while for LANS— a, 
intermittency does not vary much with a. 

Through examination of these three systems of nonlinear partial differential equations in comparison to Navier- 
Stokes, we have demonstrated that intermittency can be preserved with careful modification of the nonlinearity. 
This was seen with the LANS— a model and may be related to the conservation of small-scale circulation. Besides 
intermittency, the nonlinear terms also play a role in the energy spectrum (at both sub-filter and super-filter scales) 
and in the dissipation. These terms must model both the nonlocal interactions (to recover the intermittency) but 
also the local interactions (which are too strongly suppressed inside the "rigid bodies" of LANS— a). Finally, we have 
demonstrated that regularization modeling can be employed to reduce the computational cost while preserving the 
high-order statistics of the fiow. 

We remark that the computational gain thus far achieved by any of these regularizations is insufficient for applica- 
tions at very high Reynolds numbers, and the three subgrid stress tensors discussed here may need to be supplemented 
with an enhanced effective viscosity to be employed as LES. This is a common practice when implementing the Clark 
model (see e.g., Ref. [l^l for a study of this model with an extra Smagorinsky term). In this light, the present study 
may be useful as an analysis of the properties of the SGS tensors of the regularizations, and to pick best candidates, 
before the addition of enhanced dissipation. Studies similar to that in Ref. [4^ will also need to be done for the cases 
of Clark— a and Leray— a to quantify their computational savings before the addition of such dissipative terms. 
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